/* 
 * File:   InterpolacionNewton.cpp
 * Author: Gina Ceron
 * 
 * Created on 2 de septiembre de 2011, 09:31 AM
 */

#include "InterpolacionNewton.h"
#include "librerias/Polinomio.h"

#include <cstdlib>
#include <stdio.h>
#include <iostream>
#include <stdlib.h>
#include <string.h>
#include <stdlib.h>
#include <math.h>
#include <cstdio>

#include <set>

InterpolacionNewton::InterpolacionNewton(int numAletaorios, int tamEntrada, int nanyos, int ajuste, int ajustAleat)
{
    n_anyos=nanyos;
    tam_entrada=tamEntrada;
    tam_aleatorio=numAletaorios;
    tam_ajuste = ajuste;
    tam_ajusteAleat=ajustAleat;
    mat_DatosXY = new int* [tam_entrada];
    mat_DatosInterpol = new int* [tam_ajuste];
    mat_DiferenciasDiv1= new double* [tam_ajuste];
    mat_DiferenciasDiv2= new double* [tam_ajuste];
    mat_DifDivAleat= new double * [tam_ajusteAleat];
    mat_DatosInterAleat=new int *[tam_ajusteAleat];
    
    //Init matrizDatosInterAleat
    for(int i=0; i<tam_ajusteAleat;i++)
    {
        mat_DatosInterAleat[i]= new int [2];
        
        for(int j=0; j<2;j++)
        {
            mat_DatosInterAleat[i][j]=0;
        }
    }
    
    //Init matrizDivDivAle
    for(int i=0; i<tam_ajusteAleat;i++)
    {
        mat_DifDivAleat[i]= new double [tam_ajusteAleat];
        
        for(int j=0; j<tam_ajusteAleat;j++)
        {
            mat_DifDivAleat[i][j]=0;
        }
    }
    
    //Init matriz Datos XY
    int cols=nanyos*2;
    for(int idFilas=0; idFilas<tam_entrada; idFilas++)
    {
        mat_DatosXY[idFilas]= new int[cols];
        for(int idCols=0; idCols<cols; idCols++)
        {
            mat_DatosXY[idFilas][idCols]=0;
        }
    }
    
    //Init matrices Datos Interpolacion y diferencias Div
    for(int idFilas=0; idFilas<tam_ajuste; idFilas++)
    {
        //Init columnas de Datos a interpolar
        //

        mat_DatosInterpol[idFilas]=new int [4];
        for(int idCols=0; idCols< 4; idCols++)
        {
            mat_DatosInterpol[idFilas][idCols]=0;
        }
        //Init columnas de Diferencias Divididas
        //
        mat_DiferenciasDiv1 [idFilas]=new double [tam_ajuste];
        mat_DiferenciasDiv2 [idFilas]=new double [tam_ajuste];
        for(int idCols=0; idCols< tam_ajuste; idCols++)
        {
            mat_DiferenciasDiv1[idFilas][idCols]=0.0;
            mat_DiferenciasDiv2[idFilas][idCols]=0.0;
        }        
    }
    
    //Init array de los datos escogidos, bs y array de resultados
    //    
    arrayEscogidos = new bool [tam_entrada];
    int tamResul=tam_entrada-tam_ajuste;
    arrayResultados = new double [tamResul];
    arrayResultados2= new double [tamResul];
    for(int i=0;i<tamResul;i++)
    {
        arrayResultados2[i]=0.0;
    }
    
    for(int i=0;i<tamResul;i++)
    {
        arrayResultados[i]=0.0;
    }
    for(int i=0;i<tam_entrada;i++)
    {
        arrayEscogidos[i]=false;
    }

    
    //Init array Datos Pruebas
    mat_DatsPrueb= new int * [tam_aleatorio];
    for(int idFilas=0;idFilas<tam_aleatorio;idFilas++)
    {
        mat_DatsPrueb [idFilas]=new int [2];
        for(int idColums=0;idColums<2;idColums++)
        {
            mat_DatsPrueb[idFilas][idColums]=0;
        }
    }

}

InterpolacionNewton::~InterpolacionNewton(){

    //Delete matriz Datos XY
    //
    for(int idFilas=0; idFilas<tam_entrada; idFilas++)
    {
        delete mat_DatosXY[idFilas];
    }
    
    //Delete matrices Datos Interpolacion y diferencias Div
    //
    for(int idFilas=0; idFilas<tam_ajuste; idFilas++)
    {
        //Delete columnas de Datos a interpolar
        delete mat_DatosInterpol[idFilas];
        
        //Delete columnas de Diferencias Divididas
        delete mat_DiferenciasDiv1[idFilas];
        delete mat_DiferenciasDiv2[idFilas];        
    }
    
    
    //Delete matrices Datos Aleatorios, matrices Datos Inter Aleatorios
    //
    for(int idFilas=0; idFilas<tam_aleatorio; idFilas++)
    {
        delete mat_DatosInterAleat[idFilas];
        delete mat_DifDivAleat[idFilas];        
    }
    
    //Del matriz Pruebas
    mat_DatsPrueb= new int * [tam_aleatorio];
    for(int idFilas=0;idFilas<tam_aleatorio;idFilas++)
    {
        delete mat_DatsPrueb[idFilas];
        mat_DatsPrueb[idFilas]=NULL;
    }
    
    delete arrayResultados;
    delete arrayResultados2;
    delete array_bk;
    delete arrayEscogidos;
    delete [] mat_DatosXY;
    delete [] mat_DatosInterpol;
    delete [] mat_DiferenciasDiv1;
    delete [] mat_DiferenciasDiv2;
    delete [] mat_DifDivAleat;      
    
    n_anyos=0;
    tam_entrada=0;
    tam_aleatorio=0;
    tam_ajuste = 0;
    
    
}

bool InterpolacionNewton::enArreglo(int * array, int n, int tamano){
    for(int i=0;i<tamano;i++){
        if(array[i]==n)return true;
    }
    return false;
}

void InterpolacionNewton::inicializarDatosPrueba(){
    
    int * setX=new int[tam_aleatorio];
//    set<int>::iterator iteradorSetX=setX.begin();
    for(int i=0;i<tam_aleatorio;i++)
    {
        setX[i]=0;
    }
    srand(time(NULL)); // Instancia semilla random
    int aleatorio=0;
    bool esta=true;
    
    for(int idfils=0; idfils < tam_aleatorio; idfils++)
    {
        for (int idcols=0; idcols < 2; idcols++)
        {
            if(idcols==0)
            {
                while(esta){
                    aleatorio = (int)pow(-1, rand()%100)*rand()%200;
                    esta=enArreglo(setX,aleatorio,tam_aleatorio);
                }
                esta=true;
                setX[idfils]=aleatorio;
            }else aleatorio = (int)pow(-1, rand()%150)*rand()%800;
            
            mat_DatsPrueb[idfils][idcols]=aleatorio;
            
        }
    }
}

void InterpolacionNewton::leerArchivo(){
    FILE * ptr_archivo= fopen("datosPractica3.csv", "r+");
    if(ptr_archivo!=NULL)
    {
        int enteroTmp=0;
        char str [100];
        char str2[100];
        int cols=n_anyos*2;
        for (int i =0; i<3;i++){
        
            fscanf(ptr_archivo,"\"%s ", str);
            fscanf(ptr_archivo,"%d", &enteroTmp);
            if (i==2)fscanf(ptr_archivo,"%s\n", str2);
        }
        fscanf(ptr_archivo,"%s",str);
        for(int idFila=0;idFila<tam_entrada;idFila++)
        {
            for(int idCol=0;idCol<cols;idCol++)
            {                
                fscanf(ptr_archivo,"%d;", &enteroTmp);
                mat_DatosXY[idFila][idCol]= enteroTmp;
            }
        }
    }
    else cout<<" No existe el archivo"<<endl;
}

void InterpolacionNewton::escogerDatosInter()
{
    // Instancia semilla random
    srand(time(NULL));
    // Construyo un aleatorio para numeros entre 0 y 31
    int aleatorio;
    
    bool yaEsta = true;
    //Escojo los 20 datos escogidos a interpolar, segun el aleatorio
    for(int idfils=0; idfils < tam_ajuste; idfils++)
    {
        aleatorio=2*idfils;
        
        if(aleatorio>30){
            while(yaEsta)
            {
                aleatorio = rand()%31;
                //Chequear que no se repita aleatorios, al escoger 20 de 32 la probabilidad d repetir es alta
                yaEsta=arrayEscogidos[aleatorio];             
            }
            yaEsta=true;
            arrayEscogidos[aleatorio]=true;
        }else{
            arrayEscogidos[aleatorio]=true;
        }
        
        for(int idcols=0; idcols <4; idcols++)
        {
            mat_DatosInterpol[idfils][idcols]=mat_DatosXY[aleatorio][idcols];
            
        }
    }
}

void InterpolacionNewton::calculoDiferenciasDivididas()
{
    for(int idfils = 0; idfils<tam_ajuste; idfils++)
    {
        //Ingreso de datos x y y en matriz q guarda diferencias divididas para el 1er anho
        mat_DiferenciasDiv1[idfils][0]=(double)mat_DatosInterpol[idfils][1];
        //Ingreso de datos x y y en matriz q guarda diferencias divididas para el 2do anho
        mat_DiferenciasDiv2[idfils][0]=(double)mat_DatosInterpol[idfils][3];
    }
    
    double dividendo1, dividendo2, divisor1, divisor2=0.0;
    
    //Calculamos la tabla de diferencias divididas, para ambos años
    for(int idcols = 1; idcols<tam_ajuste; idcols++)
    {
        for(int idfils=0; idfils<tam_ajuste-idcols; idfils++)
        {
            dividendo1=(mat_DiferenciasDiv1[idfils+1][idcols-1])-(mat_DiferenciasDiv1[idfils][idcols-1]);
            dividendo2=(mat_DiferenciasDiv2[idfils+1][idcols-1])-(mat_DiferenciasDiv2[idfils][idcols-1]);
            divisor1=(double)(mat_DatosInterpol[idfils+idcols][0]-mat_DatosInterpol[idfils][0]);
            divisor2=(double)(mat_DatosInterpol[idfils+idcols][2]-mat_DatosInterpol[idfils][2]);
            mat_DiferenciasDiv1[idfils][idcols]=dividendo1/divisor1;
            mat_DiferenciasDiv2[idfils][idcols]=dividendo2/divisor2;
        }      
    }
}

void InterpolacionNewton::construccionPolinomio(){
    //Creo arreglo de polinomios por cantidad de años a calcular
    //
    array_Polinomios = new Polinomio[2];
    array_PolysError = new Polinomio[2];
    
    Polinomio polyTmp1// Tendrá almacenado el polinomio resultante de la productoria, se utilizara como el temporal para la productoria
            ,polyTmp2/* Tendrá el polinomio creado a multiplicar*/
            ,polySumaParcial/*Se utiliza para almacenar la sumatoria*/;
    long double b_diferDiv;
    long double x_j;
 
    //For para hacerlo por cada año
    //
    for(int ano=0; ano<2 ;ano++)
    {
        //Comienzo con el primer año
        //
        if(ano==0){
            b_diferDiv=(long double)mat_DiferenciasDiv1[0][0];
        }
        else b_diferDiv=(long double)mat_DiferenciasDiv2[0][0];
        
        polySumaParcial.resetPolinomio();
        polySumaParcial.crearNodoPolinomio(b_diferDiv,0);//Valor inicial, el monomio b_0
        
        //For para la sumatoria
        for(int i=1;i<tam_ajuste;i++)
        {
            polyTmp1.resetPolinomio();
            polyTmp1.crearNodoPolinomio(1.0,0);
 
            if(ano==0) b_diferDiv=mat_DiferenciasDiv1[0][i];
            else b_diferDiv=mat_DiferenciasDiv2[0][i];
            
            //For para la productoria
            //
            for(int j=0; j< i;j++)
            {
                polyTmp2.resetPolinomio();
                polyTmp2.crearNodoPolinomio(1,1);
                
                if(ano==0)x_j=(long double) mat_DatosInterpol[j][0];
                else x_j=(long double)mat_DatosInterpol[j][2];
                
                polyTmp2.crearNodoPolinomio(-x_j,0);
                polyTmp1=polyTmp1.PolinomioXPolinomio(polyTmp1,polyTmp2);
                polyTmp2.resetPolinomio();
            }
            
            polyTmp1.const_mult_poly(b_diferDiv);//multiplica la constante b por la productoria
            
            if(i<tam_ajuste-1)polySumaParcial=polySumaParcial.sumaPolinomio(polySumaParcial,polyTmp1);//Suma polinomio de productoria y polinomio parcial
            else array_PolysError[ano]=polyTmp1;        
        }
        
        array_Polinomios[ano]=polySumaParcial;
    }
    
}

void InterpolacionNewton::imprimirPolinomios(){
    for(int i=0;i<2;i++){
        cout<<"Año "<<i+1<<": ";
        array_Polinomios[i].imprimirPolinomio();
        cout<<endl;
    }
}

void InterpolacionNewton::ordenarDatosAInterpolar(){
//    Insertion
    int i, j, tmpX, tmpY;
    for(i=1;i<tam_ajuste;i++)
    {
        //Ordeno primer array
        //
        tmpX=mat_DatosInterpol[i][0];
        tmpY=mat_DatosInterpol[i][1];
        j=i-1;
        while(j >= 0)
        {
            if(mat_DatosInterpol[j][0] > tmpX)
            {
                mat_DatosInterpol[j+1][0] = mat_DatosInterpol[j][0];
                mat_DatosInterpol[j+1][1] = mat_DatosInterpol[j][1];
                j--;
            }else break;
        }
        mat_DatosInterpol[j+1][0] = tmpX;
        mat_DatosInterpol[j+1][1] = tmpY;
        
        //Ordeno segundo array
        //
        tmpX=mat_DatosInterpol[i][2];
        tmpY=mat_DatosInterpol[i][3];
        j=i-1;
        while(j >= 0)
        {
            if(mat_DatosInterpol[j][2] > tmpX)
            {
                mat_DatosInterpol[j+1][2] = mat_DatosInterpol[j][2];
                mat_DatosInterpol[j+1][3] = mat_DatosInterpol[j][3];
                j--;
            }else break;
        }
        mat_DatosInterpol[j+1][2] = tmpX;
        mat_DatosInterpol[j+1][3] = tmpY;
    }
}

void InterpolacionNewton::evaluarEnPoli(){
    
    double x;
    double y;
    int contResul=0;
    Polinomio polyTmp;
//    array_Polinomios[0].imprimirPolinomio();
//    array_Polinomios[1].imprimirPolinomio();
    for(int i=0;i<tam_entrada;i++)
    {
        if(!arrayEscogidos[i])
        {
            polyTmp=array_Polinomios[0];
            x=(long double)mat_DatosXY[i][0];
            y=polyTmp.evaluarPolinomio(x);
            arrayResultados[contResul]=y;
            polyTmp=array_Polinomios[1];
            x=(long double)mat_DatosXY[i][2];
            y=polyTmp.evaluarPolinomio(x);
            arrayResultados2[contResul]=y;
            
            contResul++;
        }
    }
}

void InterpolacionNewton::escribirResultado(){
    int contResul=0;
    FILE * output;
    output = fopen("resultadosNewton.txt","w+");
    
    int colx, coly;
    int x_entrada,y_esperado;
    double y_salida;
    
    fprintf(output,"!!!!!Punto 1!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!\n");
        
    for(int ano=0;ano<2;ano++)
    {
        fprintf(output,"!!!!!Polinomio No.%d!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!\n",ano+1);
        fprintf(output,"!!!!!X; Y; P(X);R(X); \n");
        if(ano==0)
        {
            colx=0;
            coly=1;
        }
        else{
            colx=2;
            coly=3;
        }
        for(int i=0;i<tam_entrada;i++)
        {
            x_entrada=mat_DatosXY[i][colx];
            y_esperado=mat_DatosXY[i][coly];
            if(!arrayEscogidos[i])
            {
                if(ano==0)y_salida=arrayResultados[contResul];
                else y_salida=arrayResultados2[contResul];
                contResul++;
                fprintf(output,"%d; %d; %d; %f;\n", x_entrada, y_esperado, (int)round(y_salida),(double)array_PolysError[ano].evaluarPolinomio((long double)x_entrada));

            }
            else{
                y_salida=array_Polinomios[ano].evaluarPolinomio(x_entrada);
                fprintf(output,"!!!!!%d; %d; %d; %f;\n", x_entrada, y_esperado, (int)round(y_salida),(double)array_PolysError[ano].evaluarPolinomio((long double)x_entrada));
            }
         }
        contResul=0;
    }
    fclose(output);
}

void InterpolacionNewton::escogerDatosInter2()
{
    delete arrayEscogidos;
    arrayEscogidos=NULL;
    arrayEscogidos=new bool [tam_aleatorio];
    for(int i=0;i<tam_aleatorio;i++) arrayEscogidos[i]=false;
    // Instancia semilla random
    srand(time(NULL));
    // Construyo un aleatorio para numeros entre 0 y 31
    int aleatorio;
    
    bool yaEsta = true;
    //Escojo los 20 datos escogidos a interpolar, segun el aleatorio

    for(int idfils=0; idfils < tam_ajusteAleat; idfils++)
    {
        aleatorio=2*idfils;
        
        if(aleatorio>24){
            while(yaEsta)
            {
                aleatorio = rand()%tam_aleatorio;
                //Chequear que no se repita aleatorios, al escoger 20 de 32 la probabilidad d repetir es alta
                yaEsta=arrayEscogidos[aleatorio];             
            }
            yaEsta=true;
            arrayEscogidos[aleatorio]=true;
        }else{
            arrayEscogidos[aleatorio]=true;
        }
        for(int idcols=0; idcols <2; idcols++)
        {
            mat_DatosInterAleat[idfils][idcols]=mat_DatsPrueb[aleatorio][idcols];
        }
    }

}

void InterpolacionNewton::ordenarDatosInterPAleato(){
    
//    Insertion
    int i, j, tmpX, tmpY;
    for(i=1;i<tam_ajusteAleat;i++)
    {
        //Ordeno primer array
        //
        tmpX=mat_DatosInterAleat[i][0];
        tmpY=mat_DatosInterAleat[i][1];
        j=i-1;
        while(j >= 0)
        {
            if(mat_DatosInterAleat[j][0] > tmpX)
            {
                mat_DatosInterAleat[j+1][0] = mat_DatosInterAleat[j][0];
                mat_DatosInterAleat[j+1][1] = mat_DatosInterAleat[j][1];
                j--;
            }else break;
        }
        mat_DatosInterAleat[j+1][0] = tmpX;
        mat_DatosInterAleat[j+1][1] = tmpY;
    }    
}

void InterpolacionNewton::calculoDiferenciasDivididasSegundo()
{
    for(int idfils = 0; idfils<tam_ajusteAleat; idfils++)
    {
        //Ingreso de datos x y y en matriz q guarda diferencias divididas para el 1er anho
        mat_DifDivAleat[idfils][0]=(double)mat_DatosInterAleat[idfils][1];
    }
    
    double dividendo, divisor;
    
    //Calculamos la tabla de diferencias divididas, para ambos años
    for(int idcols = 1; idcols<tam_ajusteAleat; idcols++)
    {
        for(int idfils=0; idfils<tam_ajusteAleat-idcols; idfils++)
        {
            dividendo=(mat_DifDivAleat[idfils+1][idcols-1])-(mat_DifDivAleat[idfils][idcols-1]);
            divisor=(double)(mat_DatosInterAleat[idfils+idcols][0]-mat_DatosInterAleat[idfils][0]);
            mat_DifDivAleat[idfils][idcols]=dividendo/divisor;
        }      
    }
}

void InterpolacionNewton::construccionPolinomio2(){
    //Creo arreglo de polinomios por cantidad de años a calcular
    //
    Polinomio polyTmp1// Tendrá almacenado el polinomio resultante de la productoria, se utilizara como el temporal para la productoria
            ,polyTmp2/* Tendrá el polinomio creado a multiplicar*/
            ,polySumaParcial/*Se utiliza para almacenar la sumatoria*/;
    long double b_diferDiv;
    long double x_j;
 
    //Comienzo con el primer año
    //
    b_diferDiv=(long double)mat_DifDivAleat[0][0];
    
    polySumaParcial.resetPolinomio();
    polySumaParcial.crearNodoPolinomio(b_diferDiv,0);//Valor inicial, el monomio b_0

    //For para la sumatoria
    for(int i=1;i<tam_ajusteAleat;i++)
    {
        polyTmp1.resetPolinomio();
        polyTmp1.crearNodoPolinomio(1.0,0);

        b_diferDiv=mat_DifDivAleat[0][i];

        //For para la productoria
        //
        for(int j=0; j< i;j++)
        {
            polyTmp2.resetPolinomio();
            polyTmp2.crearNodoPolinomio(1,1);

            x_j=(long double)mat_DatosInterAleat[j][0];

            polyTmp2.crearNodoPolinomio(-x_j,0);
            polyTmp1=polyTmp1.PolinomioXPolinomio(polyTmp1,polyTmp2);
        }

        polyTmp1.const_mult_poly(b_diferDiv);//multiplica la constante b por la productoria

        if(i<tam_ajuste-1)polySumaParcial=polySumaParcial.sumaPolinomio(polySumaParcial,polyTmp1);//Suma polinomio de productoria y polinomio parcial
        else polyError2=polyTmp1;        
    }

    poly2=polySumaParcial;
}

void InterpolacionNewton::evaluarEnPoliSegundo(){
    delete arrayResultados2;
    arrayResultados=NULL;
    arrayResultados2=new double [tam_ajusteAleat];
    double x;
    double y;
    int contResul=0;
//    array_Polinomios[0].imprimirPolinomio();
//    array_Polinomios[1].imprimirPolinomio();
    for(int i=0;i<tam_aleatorio;i++)
    {
        if(!arrayEscogidos[i])
        {
            x=(long double)mat_DatsPrueb[i][0];
            y=poly2.evaluarPolinomio(x);
            arrayResultados2[contResul]=y;
            contResul++;
        }
    }
}

void InterpolacionNewton::escribirPuntoSegundo(){
    int contResul=0;
    FILE * output;
    output = fopen("resultadosNewton2.txt","w+");
    
    int colx=0, coly=1;
    int x_entrada,y_esperado;
    double y_salida;
    
    fprintf(output,"!!!!!Punto 2!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!\n");
        
    
    fprintf(output,"!!!!!Polinomio No.%d!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!\n");
    fprintf(output,"!!!!!X; Y; P(X);R(X); \n");
    for(int i=0;i<tam_aleatorio;i++)
    {
        x_entrada=mat_DatsPrueb[i][colx];
        y_esperado=mat_DatsPrueb[i][coly];
        if(!arrayEscogidos[i])
        {
            y_salida=arrayResultados2[contResul];
            contResul++;
            fprintf(output,"%d; %d; %d; %f;\n", x_entrada, y_esperado, (int)round(y_salida),(double)polyError2.evaluarPolinomio((long double)x_entrada));

        }
        else{
            y_salida=poly2.evaluarPolinomio(x_entrada);
            fprintf(output,"!!!!!%d; %d; %d; %f;\n", x_entrada, y_esperado, (int)round(y_salida),(double)polyError2.evaluarPolinomio((long double)x_entrada));
        }
     }
    
    fclose(output);    
}